function conditionNumber = calculateConditionNumber3DByPowerMethod( H, ErrTol, MaxIts, verbose,n,flag,precision)
    if flag == 0
        L=3*n-2;
        allZeros=zeros(L,L,L,precision);
        Hhat=fftn(H.filter,[L,L,L]);
    end
    g=zeros(n,n,n,precision);
	W=randn(n,n,n,precision);
	W=W/norm(W(:));
    v=W;
    k=100;    
    for i = 1:k
        if flag == 0
            v=applyHFilter3D(v,Hhat);
        else
            v=optimizedAdjointPPFT3D(optimizedPPFT3D(v,precision),precision);
        end
        d=norm(v(:));
        rho=sqrt(d);
        v=v/d;
    end
    normH=rho;
    v=W;
    for i = 1:k
        if flag == 0
            [v,~,~,~] = ConjugateGradient3D(v,Hhat,n,allZeros,ErrTol,MaxIts,g,verbose);
        else
            %note H is a preconditioner, not a filter kernel
            [v,~,~,~] = ConjugateGradientFeval('diagonalPrecon3DApplication',v,{precision,H,n},ErrTol,MaxIts,g,verbose);
        end
        d=norm(v(:));
        rho=sqrt(d);
        v=v/d;
    end
    normHinv=rho;
    conditionNumber=(normH*normHinv)^2;
end
